function alpha = solveAlpha(alpha0,lambda,sigma,tspan,dt)
% alpha = solveAlpha(alpha0,lambda,sigma,t,dt) integrate a stochastic ODE
% and give the value of alpha(t) at time points defined in tspan.
% alpha0 is a scalar
% lambda, sigma is scalar
% tspan is a (P+1)-by-1 vector
% alpha is a P-by-1 vector
alpha = alpha0;
t = tspan(1);
sqrt_dt = sqrt(dt);
while t<tspan(2)
    alpha = alpha+lambda*alpha*dt+sigma*randn(1,1)*sqrt_dt;
    t = t+dt;
end